scr <- read.dta("~/Dropbox/diffusion vs. osmosis/3 Data and Code/green taxes and subsidies/OECD subsidies.dta")
library(systemfit)
variable<-c("obs_in", "ccode", "year","rgdppc","unemp", "unempsq", "i_itax_per_cap", "actual_economic_flows", "rile_const_w","env_const_w",
"env_const_w_sq","green_party","en_prod_rgdp","country_2",
"country_3",
"country_4", "country_5",
"country_6","country_7", "country_8", "country_9",  "country_10", "country_11"
,"country_12","country_13","country_14","country_15","country_16","country_17","country_18","country_19",
"country_20",
"country_21","country_22","country_23","country_24","country_25","country_26","country_27"
)
# names of all the variables in the regression: including the DV;
# reduce the original data set to a parsimomious one having only variables included in the regression;
dat.p<-subset(scr, select=variable)
dat.p$ccode<-as.factor(as.character(dat.p$ccode))
# indicator vector for those rows having missing values:
NO<-NULL
for (i in 1: dim(dat.p)[1]) {if (is.na(rowSums(dat.p[,-1]))[i]){NO<-c(NO, i)}}
dat.p<-dat.p[-NO,]
################################################################################
getwd() # Adjust directory
# Load packages:
library(foreign)
library(calibrate)
require(stats)
library(Hmisc)
library(MCMCpack)
library(ggplot2)
library(austin)
#install cshapes package
install.packages("cshapes", dependencies=TRUE)
#call/load the cshape package (and others)
library(cshapes)
library(sp)
library(spdep)
install.packages("spdep")
library(spdep)
#Retrieve the dataset (all data)
cshp.data <- cshp()
install.packages("maps")
?maps
??maps
install.packages("knitr")
install.packages("knitcitations")
setwd("/Users/genovesefederica/Dropbox/cc_encycl/encyclicals_model fit output")
library(foreign)
library(tm)
library(wordcloud)
library(SnowballC)
library(austin)
library(lda)
library(reshape)
library(reshape2)
library(modeltools)
library(ggplot2)
library(topicmodels)
library(dplyr)
library(grid)
library(Hmisc)
library(rJava)
#library(openNLP)
library(base)
library(lda)
library(plyr)
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(2900,1600),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="blue",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(2900,1600),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4, "5"),las=01,cex.axis=.85)
axis(2, at=c(1600,1900,2200,2500,2900), lab=c("1600",  "1900",  "2200","2500","2900"), las = 0,cex.axis=.85)
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(1600,2900),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="blue",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(2900,1600),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4, "5"),las=01,cex.axis=.85)
axis(2, at=c(1600,1900,2200,2500,2900), lab=c("1600",  "1900",  "2200","2500","2900"), las = 0,cex.axis=.85)
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(1600,2900),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="blue",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(1600,2900),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4, "5"),las=01,cex.axis=.85)
axis(2, at=c(1600,1900,2200,2500,2900), lab=c("1600",  "1900",  "2200","2500","2900"), las = 0,cex.axis=.85)
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(1600,2900),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="blue",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(1600,2900),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4", "5"),las=01,cex.axis=.85)
axis(2, at=c(1600,1900,2200,2500,2900), lab=c("1600",  "1900",  "2200","2500","2900"), las = 0,cex.axis=.85)
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(1600,2900),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="blue",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(1600,2900),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4", "5"),las=01,cex.axis=.85)
axis(2, at=c(1600,1900,2200,2500,2900), lab=c("1600",  "1900",  "2200","2500","2900"), las = 0,cex.axis=.85)
#Next plot loglike measure
combinell<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,4], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(-25000000,-23000000),main="")
par(new=T)
combinell<-cbind(combinell,forplot[,4])
}
plot(rowMeans(combinell), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Log-Likelihood",axes=FALSE,ylim=c(-25000000,-23000000),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4, "5"),las=01,cex.axis=.85)
axis(2, at=c(-25000000, -24500000,-24000000, -23500000, -23000000), lab=c( "-25000000", "-24500000","-24000000",  "-23500000","-23000000"), las = 0,cex.axis=.85)
#Next plot loglike measure
combinell<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,4], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(-25000000,-23000000),main="")
par(new=T)
combinell<-cbind(combinell,forplot[,4])
}
plot(rowMeans(combinell), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Log-Likelihood",axes=FALSE,ylim=c(-25000000,-23000000),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4", "5"),las=01,cex.axis=.85)
axis(2, at=c(-25000000, -24500000,-24000000, -23500000, -23000000), lab=c( "-25000000", "-24500000","-24000000",  "-23500000","-23000000"), las = 0,cex.axis=.85)
#Next plot loglike measure
combinell<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,4], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(-2500000,-2300000),main="")
par(new=T)
combinell<-cbind(combinell,forplot[,4])
}
plot(rowMeans(combinell), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Log-Likelihood",axes=FALSE,ylim=c(-2500000,-2300000),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4", "5"),las=01,cex.axis=.85)
axis(2, at=c(-2500000, -2450000,-2400000, -2350000, -2300000), lab=c( "-2500000", "-2450000","-2400000",  "-2350000","-2300000"), las = 0,cex.axis=.85)
combinell<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,4], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(-2500000,-2300000),main="")
par(new=T)
combinell<-cbind(combinell,forplot[,4])
}
plot(rowMeans(combinell), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Log-Likelihood",axes=FALSE,ylim=c(-2500000,-2300000),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4", "5"),las=01,cex.axis=.85)
axis(2, at=c(-2500000, -2400000, -2300000), lab=c( "-2500000","-2400000", "-2300000"), las = 0,cex.axis=.85)
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(1600,2900),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="blue",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(1600,2900),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4", "5"),las=01,cex.axis=.85)
axis(2, at=c(1600,1900,2200,2500,2900), lab=c("1600",  "1900",  "2200","2500","2900"), las = 0,cex.axis=.85)
#Next plot loglike measure
combinell<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,4], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(-2500000,-2300000),main="")
par(new=T)
combinell<-cbind(combinell,forplot[,4])
}
plot(rowMeans(combinell), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Log-Likelihood",axes=FALSE,ylim=c(-2500000,-2300000),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4", "5"),las=01,cex.axis=.85)
axis(2, at=c(-2500000, -2400000, -2300000), lab=c( "-2500000","-2400000", "-2300000"), las = 0,cex.axis=.85)
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfitfull",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(1600,2900),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(1600,2900),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4", "5"),las=01,cex.axis=.85)
axis(2, at=c(1600,1900,2200,2500,2900), lab=c("1600",  "1900",  "2200","2500","2900"), las = 0,cex.axis=.85)
#Sys.setenv(NOAWT= "true") # Needed for clean up of texts
library(foreign)
library(tm)
library(wordcloud)
library(SnowballC)
library(austin)
library(lda)
library(reshape)
library(reshape2)
library(modeltools)
library(ggplot2)
library(topicmodels)
library(dplyr)
library(grid)
library(Hmisc)
library(rJava)
#library(openNLP)
library(base)
library(lda)
library(plyr)
getwd() #Fix based on your directory
setwd("/Users/genovesefederica/Dropbox/cc_encycl/encyclicals_it")
cname <- file.path(".", "encyclicals_it_re")
cname
dir(cname)
#convert text list  to corpus
docs <- Corpus(DirSource(cname))
class(docs)
class(docs[[1]])
summary(docs)
inspect(docs[1])
encycl <- rownames(summary(docs))
encycl
getTransformations()
# replace "/" used sometimes to separate alternative words
for (j in seq(docs))
{
docs[[j]] <- gsub("/", " ", docs[[j]])
docs[[j]] <- gsub("@", " ", docs[[j]])
}
# Preparation of documents (conversion to lower case, removal of punctuation, yada yada)
docs <- tm_map(docs, removeNumbers)
docs <- tm_map(docs, tolower)
docs <- tm_map(docs, PlainTextDocument)
docs <- tm_map(docs, removeNumbers)
docs <- tm_map(docs, removePunctuation)
docs <- tm_map(docs, stripWhitespace)
docs2 <- tm_map(docs, stemDocument)
docs <- tm_map(docs2, removeWords, c(stopwords('italian'),"che", "ma", "per", "anch", "ogni", "modo", "piu", "perch", "cos","alcuni", "sempr", "quando", "cose", "part", "senza", "esser","aver","cio","puo","possono","molto","tutt","fatto", "solo","tale","stess","pid","qual","gesd","deg","deve","infatti","dell","essa","stesso","ag","ci","gia","delluomo","eg","cf","pu","d"))
#docs <- docs
inspect(docs[1])
#create a document term matrix object for input into the LDA model
dtm <- DocumentTermMatrix(docs)
dtm
#Keep some important elements of this matrix for post-estimation routines
words<-dtm[[6]]
document_index<-dtm[[4]]
# Summary stats
rows <- rowSums(as.matrix(dtm))
rows
smean.sd(rows)
smedian.hilow(rows, conf.int=.05)
# Check out the raw frequencies
m <- as.matrix(dtm)
v <- sort(colSums(m), decreasing=TRUE)
head(v, 70)
words <- names(v)
d <- data.frame(word=words, freq=v)
#wordcloud(d$word, d$freq, min.freq=500)
wordcloud(d$word, d$freq, min.freq=500, random.order=FALSE, rot.per=0)
#df <- d[d$freq >= 500, ]
rownames(m) <- encycl
newnames <- c("1959 ad petri ","1959 sacerdotii ","1959 grata recordatio ","1959 princeps ","1961 mater ","1961 aeterna dei ","1962 paenitentiam ","1963 pacem ","1964 ecclesiam ","1965 mense maio ","1965 mysterium ","1966 christi matri ","1967 populorum progressio ","1967 sacerdotalis caelibatus ","1968 humanae vitae ","1979 redemptor hominis ","1980 dives in misericordia ","1981 laborem exercens ","1985 slavorum apostoli ","1986 dominum et vivificantem ","1987 redemptoris mater ","1987 sollicitudo rei socialis ","1990 redemptoris missio ","1991 centesimus annus ","1993 veritatis splendor ","1995 evangelium vitae ","1995 ut unum sint ","1998 fides et ratio ","2003 eccl de euch ","2005 deus caritas est ","2007 spe salvi ","2009 caritas in veritate ","2013 enciclica lumen fidei it ","2015 laudato si ")
rownames(m) <- newnames
m2 <- t(m)
new_m <- wfm(m2)
dim(new_m)
docs(new_m)
wfm.to.lda <- function(wfm){
m <- as.worddoc(wfm)
v <- words(m)
d <- list()
for (i in 1:length(docs(m))){
nzero <- which(m[,i]>0)
d[[i]] <- t(matrix(as.integer(c(nzero-1, m[nzero, i])), ncol=2))
}
t <- docs(m)
corpus <- list(vocab=v, documents=d, titles=t)
return(corpus)
}
wfm.to.bmr <- function (y, wfm, filename, line.sep='\n') {
x <- as.worddoc(wfm)
lines <- c()
marg <- wordmargin(wfm)
for (i in 1:ncol(x)) {
nonz.vals <- c()
if (marg == 1){
## docs are columns
nonz <- which(x[, i] != 0)
if (length(nonz)>0)
nonz.vals <- as.numeric(x[nonz, i])
} else {
## docs are rows
nonz <- which(x[i, ] != 0)
if (length(nonz)>0)
nonz.vals <- as.numeric(x[i, nonz])
}
if (length(nonz.vals)>0){
## construct a line
## 0 if we don't know the class (ignored in training set)
yy <- ifelse(!is.null(y), y[i], 0)
ldat <- paste(nonz, ":", nonz.vals, collapse = " ", sep = "")
lines <- c(lines, paste(yy, ldat))
}
}
writeLines(lines, filename, sep=line.sep)
}
# 2 topics
m.all <- wfm(m, word.margin=2)
m.lda <- wfm.to.lda(m.all)
mod <- lda.collapsed.gibbs.sampler(m.lda$documents, 2,
m.lda$vocab, 50, 1, 1, compute.log.likelihood=TRUE)
# Show the topic model results
top.topic.words(mod$topics, 10, by.score=TRUE)
top.words <- top.topic.words(mod$topics, 10, by.score=TRUE)
top.words
topics <-  c(1, 2, 3, 4, 5)
perplex<-NULL
loglike<-NULL
for (k in topics) {
lda <- LDA(m, control=list(seed = 123, burnin = 500, thin = 100, iter = 1000), k = k, compute.log.likelihood = TRUE)
perplex<-rbind(perplex,perplexity<-perplexity(lda, m))
loglike<-rbind(loglike,logLik(lda))
}
topics <-  c(2, 3, 4, 5, 6)
for (k in topics) {
lda <- LDA(m, control=list(seed = Seed, burnin = 500, thin = 100, iter = 1000), k = k, compute.log.likelihood = TRUE)
perplex<-rbind(perplex,perplexity<-perplexity(lda, m))
loglike<-rbind(loglike,logLik(lda))
}
for (k in topics) {
lda <- LDA(m, control=list(seed = 123, burnin = 500, thin = 100, iter = 1000), k = k, compute.log.likelihood = TRUE)
perplex<-rbind(perplex,perplexity<-perplexity(lda, m))
loglike<-rbind(loglike,logLik(lda))
}
for (k in topics) {
lda <- LDA(m, control=list( burnin = 500, thin = 100, iter = 1000), k = k, compute.log.likelihood = TRUE)
perplex<-rbind(perplex,perplexity<-perplexity(lda, m))
loglike<-rbind(loglike,logLik(lda))
}
for (k in topics) {
lda <- LDA(m, control=list(seed = 123), k = k, compute.log.likelihood = TRUE)
perplex<-rbind(perplex,perplexity<-perplexity(lda, m))
loglike<-rbind(loglike,logLik(lda))
}
forplot<-cbind(topics,perplex,loglike)
#save combined results for subsequent plotting
forplot<-as.data.frame(forplot)
write.csv(forplot,paste("modelfit",i,".csv"))
getwd()
setwd("/Users/genovesefederica/Dropbox/cc_encycl/encyclicals_model fit output")
for(i in 1:10){
#start with null perplexity and loglike objects
perplex<-NULL
loglike<-NULL
#Generate range of topic-numbers to be evaluated
topics <-  c(2, 3, 4, 5, 6)
#set the random number seed for this iteration
Seed<-i
for (k in topics) {
lda <- LDA(m, control=list(seed = Seed, burnin = 500, thin = 100, iter = 1000), k = k, compute.log.likelihood = TRUE)
perplex<-rbind(perplex,perplexity<-perplexity(lda, m))
loglike<-rbind(loglike,logLik(lda))
}
forplot<-cbind(topics,perplex,loglike)
#save combined results for subsequent plotting
forplot<-as.data.frame(forplot)
write.csv(forplot,paste("modelfit",i,".csv"))
} #end loop
for(i in 1:10){
#start with null perplexity and loglike objects
perplex<-NULL
loglike<-NULL
#Generate range of topic-numbers to be evaluated
topics <-  c(2, 3, 4, 5, 6)
#set the random number seed for this iteration
Seed<-i
for (k in topics) {
lda <- LDA(m, control=list(seed = Seed), k = k, compute.log.likelihood = TRUE)
perplex<-rbind(perplex,perplexity<-perplexity(lda, m))
loglike<-rbind(loglike,logLik(lda))
}
forplot<-cbind(topics,perplex,loglike)
#save combined results for subsequent plotting
forplot<-as.data.frame(forplot)
write.csv(forplot,paste("modelfit",i,".csv"))
} #end loop
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfit",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(3200,4500),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(3200,4500),main="")
axis(1, at=1:5, lab=c("2", "3",  "4",  "5", "6"),las=01,cex.axis=.85)
axis(2, at=c(3200,3850,4500), lab=c("3200",    "3850","4500"), las = 0,cex.axis=.85)
for(i in 1:10){
#start with null perplexity and loglike objects
perplex<-NULL
loglike<-NULL
#Generate range of topic-numbers to be evaluated
topics <-  c(2, 3, 4, 5, 6)
#set the random number seed for this iteration
Seed<-i
for (k in topics) {
lda <- LDA(m, control=list(seed = Seed), k = k, compute.log.likelihood = TRUE)
perplex<-rbind(perplex,perplexity<-perplexity(lda, m))
loglike<-rbind(loglike,logLik(lda))
}
forplot<-cbind(topics,perplex,loglike)
#save combined results for subsequent plotting
forplot<-as.data.frame(forplot)
write.csv(forplot,paste("modelfit",i,".csv"))
} #end loop
#Next plot loglike measure
combinell<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfit",i,".csv"))
plot(forplot[,4], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(-2500000,-2300000),main="")
par(new=T)
combinell<-cbind(combinell,forplot[,4])
}
plot(rowMeans(combinell), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Log-Likelihood",axes=FALSE,ylim=c(-2500000,-2300000),main="")
axis(1, at=1:5, lab=c("1", "2",  "3",  "4", "5"),las=01,cex.axis=.85)
axis(2, at=c(-2500000, -2400000, -2300000), lab=c( "-2500000","-2400000", "-2300000"), las = 0,cex.axis=.85)
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfit",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(2700,5300),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(2700,5300),main="")
axis(1, at=1:5, lab=c("2", "3",  "4",  "5", "6"),las=01,cex.axis=.85)
axis(2, at=c(2700,4000,5300), lab=c("2700",    "4000","5300"), las = 0,cex.axis=.85)
#Next plot loglike measure
combinell<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfit",i,".csv"))
plot(forplot[,4], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(-2700000,-2500000),main="")
par(new=T)
combinell<-cbind(combinell,forplot[,4])
}
plot(rowMeans(combinell), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Log-Likelihood",axes=FALSE,ylim=c(-2700000,-2500000),main="")
axis(1, at=1:5, lab=c("2", "3",  "4",  "5", "6"),las=01,cex.axis=.85)
axis(2, at=c(-2700000, -2600000, -2500000), lab=c( "-2700000","-2600000", "-2500000"), las = 0,cex.axis=.85)
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfit",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(2700,5300),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(2700,5300),main="")
axis(1, at=1:5, lab=c("2", "3",  "4",  "5", "6"),las=01,cex.axis=.85)
axis(2, at=c(2700,4000,5300), lab=c("2700",    "4000","5300"), las = 0,cex.axis=.85)
#Next plot loglike measure
combinell<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfit",i,".csv"))
plot(forplot[,4], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(-2700000,-2500000),main="")
par(new=T)
combinell<-cbind(combinell,forplot[,4])
}
plot(rowMeans(combinell), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Log-Likelihood",axes=FALSE,ylim=c(-2700000,-2500000),main="")
axis(1, at=1:5, lab=c("2", "3",  "4",  "5", "6"),las=01,cex.axis=.85)
axis(2, at=c(-2700000, -2600000, -2500000), lab=c( "-2700000","-2600000", "-2500000"), las = 0,cex.axis=.85)
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfit",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(2700,5300),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(2700,5300),main="")
axis(1, at=1:5, lab=c("2", "3",  "4",  "5", "6"),las=01,cex.axis=.85)
axis(2, at=c(2700,4000,5300), lab=c("2700",    "4000","5300"), las = 0,cex.axis=.85)
#First plot perplexity measure
combineperplex<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfit",i,".csv"))
plot(forplot[,3], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(2700,5300),main="")
par(new=T)
combineperplex<-cbind(combineperplex,forplot[,3])
}
plot(rowMeans(combineperplex), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Perplexity",axes=FALSE,ylim=c(2700,5300),main="")
axis(1, at=1:5, lab=c("2", "3",  "4",  "5", "6"),las=01,cex.axis=.85)
axis(2, at=c(2700,4000,5300), lab=c("2700",    "4000","5300"), las = 0,cex.axis=.85)
#pdf('ModelFitPerplexity.pdf')
pdf('ModelFitPerplexity.pdf')
dev.off()
#Next plot loglike measure
combinell<-NULL
for(i in 1:10){
forplot<-read.csv(paste("modelfit",i,".csv"))
plot(forplot[,4], type="o", col="lightgrey",xlab="",ylab="",axes=FALSE,ylim=c(-2700000,-2500000),main="")
par(new=T)
combinell<-cbind(combinell,forplot[,4])
}
plot(rowMeans(combinell), type="o", col="red",lwd=2.5,xlab="Topics",ylab="Log-Likelihood",axes=FALSE,ylim=c(-2700000,-2500000),main="")
axis(1, at=1:5, lab=c("2", "3",  "4",  "5", "6"),las=01,cex.axis=.85)
axis(2, at=c(-2700000, -2600000, -2500000), lab=c( "-2700000","-2600000", "-2500000"), las = 0,cex.axis=.85)
